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ABSTRACT 

We derive analytical formulae for the Minkowski Functions of the cosmic microwave background (CMB) 
and large-scale structure (LSS) from primordial non-Gaussianity. These formulae enable us to estimate a 
non-linear coupling parameter, /nl, directly from the CMB and LSS data without relying on numerical 
simulations of non-Gaussian primordial fluctuations. One can use these formulae to estimate statistical 
errors on /nl from Gaussian realizations, which are much faster to generate than non-Gaussian ones, 
fully taking into account the cosmic/sampling variance, beam smearing, survey mask, etc. We show 
that the CMB data from the Wilkinson Microwave Anisotropy Probe should be sensitive to |/nl| — 40 
at the 68% confidence level. The Planck data should be sensitive to |/nl| — 20. As for the LSS data, 
the late-time non-Gaussianity arising from gravitational instability and galaxy biasing makes it more 
challenging to detect primordial non-Gaussianity at low redshifts. The late-time effects obscure the 
primordial signals at small spatial scales. High-redshift galaxy surveys at z > 2 covering ^ 10 Gpc^ 
volume would be required for the LSS data to detect |/nl| — 100. Minkowski Functionals are nicely 
complementary to the bispectrum because the Minkowski Functionals are defined in real space and the 
bispectrum is defined in Fourier space. This property makes the Minkowski Functionals a useful tool in 
the presence of real-world issues such as anisotropic noise, foreground and survey masks. Our formalism 
can be extended to scale-dependent /nl easily. 

Subject headings: cosmic microwave background - large-scale structure of universe - methods: analytical 



1. INTRODUCTION 

Recent observations of cosmological fluctuations from 
the Cosmic Microwave Background (CMB) and Large 
Scale Structure (LSS) strongly support basic predictions of 
inflationary scenarios: primordial fluctuations are nearly 
scale-invariant (Spergel et al. 2003; Tegmark et al. 2004; 
Seljak et al. 2005; Spergel et al. 2006), adiabatic (Peiris et 
al. 2003; Bucher et al. 2004; Bean et al. 2006), and Gaus- 
sian (Komatsu et al. 2002, 2003; CremineUi et al. 2006; 
Spergel et al. 2006, and references therein). In order to 
discriminate between more-than-100 candidate inflation- 
ary models, however, one needs to look for deviations from 
the scale-invariance, adiabaticity as well as Gaussianity, 
for which different inflationary models make specific pre- 
dictions. Inflationary models based upon a slowly rolling 
single-field scalar field generally predict very small devia- 
tions from Gaussianity; however, the post-inflationary evo- 
lution of non-linear metric perturbations inevitably gener- 
ates ubiquitous non-Gaussian fluctuations. On the other 
hand, a broad class of inflationary models based upon dif- 
ferent assumptions about the nature of scalar fleld(s) can 
generate significant primordial non-Gaussianity (Lyth et 
al. 2003; Dvah et al. 2004; Arkami-Hamed et al. 2004; Al- 
ishahiha et al. 2004; Bartolo et al. 2004). Therefore, Gaus- 
sianity of the primordial fluctuations offers a direct test of 
inflationary models. 

It is customary to adopt the following simple form of 
primordial non-Gaussianity in Bardeen's curvature pertur- 
bations during the matter era (e.g., Komatsu & Spergel 
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2001): 



<f-0 + /NL(0'-(0')), 



(1) 



where </> is an auxiliary random-Gaussian field and /nl 
characterizes the amplitude of a quadratic correction to 
the curvature perturbations. Note that $ is related to 
the primordial comoving curvature perturbations gener- 
ated during inflation, 7?., by $ = (3/5)7?.. While this 
quadratic form is motivated by inflationary models based 
upon a single slowly-rolling scalar field, the actual pre- 
dictions usually include momentum dependence in /nl- 
(That is to say, /nl is not a constant.) Therefore, when 
precision is required, one should use the actual formula 
given by specific processes, either from primordial non- 
Gaussianity from inflation or the post-inflationary evolu- 
tion of non-linear perturbations, to calculate a more ac- 
curate form of statistical quantities such as the angular 
bispectrum of CMB (Babich et al. 2004; Liguori et al. 
2006). Nevertheless, a constant /nl is still a useful param- 
eterization of non-Gaussianity which enables us to obtain 
simple analytical formulae for the statistical quantities to 
compare with observations. The use of a constant /nl is 
also justified by the fact that the current observations are 
not sensitive enough to detect momentum-dependence of 
/nl- Therefore, we adopt the constant /nl for our analysis 
throughout this paper. Note that it is actually straightfor- 
ward to extend our formalism to a momentum-dependent 
/nl- 

So far, analytical formulae for the statistical quantities 
of the CMB from primordial non-Gaussianity are known 
only for the angular bispectrum (Komatsu & Spergel 
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2001; Babich & Zaldarriaga 2004; Liguori et al. 2006) 
and trispectrum (Okainoto & Hu 2002; Kogo & Komatsu 
2006). The analytical formulae are extremely valuable 
especially when one tries to measure non-Gaussian sig- 
nals from the data. Fast, nearly optimal estimators for 
/nl have been derived on the basis of these analytical 
formulae (Komatsu et al. 2005; Creminelli et al. 2006), 
and have been successfully applied to the CMB data from 
the Wilkinson Microwave Anisotropy Probe (WMAP): the 
current constraint on /nl from the angular bispcctrum is 
—54 to 114 at the 95% confidence level (Komatsu et al. 
2003; Spergel et al. 2006). (See CremineUi et al. 2006, for 
an alternative parameterization of /nl-) As for the LSS, 
the analytical formula is known only for the 3-d bispec- 
trum (Vcrdc ct al. 2000; Scoccimarro ot al. 2004). The 
LSS bispectrum contains not only the primordial non- 
Gaussianity, but also the late-time non-Gaussianity from 
gravitational instability and galaxy biasing, which poten- 
tially obscure the primordial signatures. 

In this paper, we derive analytical formulae for an- 
other statistical tool, namely the Minkowski Functionals 
(MPs), which describe morphological properties of fluctu- 
ating fields (Mcckc ct al. 1994; Schmalzing & Buchcrt 1997; 
Schmalzing & Gorski 1998; Winitzki & Kosowsky 1998). 
In d-dimensional space {d = 2 for CMB and d — 3 for LSS), 
d + 1 MPs are defined, as listed in Table 1. The "Euler 
characteristic" measures topology of the fields, and is es- 
sentially given by the number of hot spots minus the num- 
ber of cold spots when d = 2. This quantity is sometimes 
called the "genus statistics" , which was independently re- 
discovered by Gott et al. (1986) in search of a topological 
measure of non-Gaussianity in the cosmic density fields. 
(The Elder characteristic and genus are different only by 
a numerical coefficient, —1/2.) 

Why study MPs? Since different statistical methods 
are sensitive to different aspects of non-Gaussianity, one 
should study as many statistical methods as possible. 
Most importantly, the MPs and bispectrum are very differ- 
ent in that MPs are defined in real space, whereas the bis- 
pectrum is defined in Fourier (or harmonic) space. There- 
fore, these statistical methods are nicely complementary 
to each other. Previously there are several attempts to 
give constraints on the primordial non-Gaussianity using 
MPs (e.g, Novikov et al. 2000). Although we shall show in 
this paper that the MPs do not contain information more 
than the bispectrum in the limit that non-Gaussianity is 
weak, the complementarity is still powerful in the pres- 
ence of complicated real-world issues such as inhomoge- 
neous noise, survey mask, foreground contamination, etc. 
The MPs have also been used to constrain /nl- Komatsu 
et al. (2003) and Spergel et al. (2006) have used numeri- 
cal simulations of non-Gaussian CMB sky maps to calcu- 
late the predicted form of MPs as a function of /nl, and 
compared the predictions with the WMAP data to con- 
strain /nl, obtaining similar constraints to the bispectrum 
ones. This method (calculating the form of MPs from non- 
Gaussian simulations) is, however, a painstaking process: 
it takes about three hours to simulate one non-Gaussian 
map on one processor of SGI Origin 300. When cosmo- 
logical parameters are varied, one needs to re-simulate a 
whole batch of non-Gaussian maps from the beginning — 
this is a highly inefficient approach. Once we have the an- 



alytical formula for the MPs as a function of /nl, however, 
we no longer need to simulate non-Gaussian maps, greatly 
speeding up the measurement of /nl from the data. 

We use the perturbative formula for MPs originally de- 
rived by Matsubara (1994, 2003): assuming that non- 
Gaussianity is weak, which has been justified by the cur- 
rent constraints on /nl, we consider the lowest-order cor- 
rections to the MPs using the multi-dimensional Edge- 
worth expansion around a Gaussian distribution function. 

The organization of paper is as follows; In § 2 we review 
the generic perturbative formula for the Minkowski Punc- 
tionals. In § 3 we derive the analytical formula for MPs of 
the CMB from primordial non- Gaussian fluctuations pa- 
rameterized by /nl- We also estimate projected statistical 
errors on /nl expected from the WMAP data from multi- 
year observations as well as from the Planck data. In § 4 
we derive the analytical formula for MPs of the LSS from 
primordial non-Gaussianity, non-linear gravitational evo- 
lution, and galaxy biasing in a perturbative manner. § 5 is 
devoted to summary and conclusions. In Appendix A we 
outline our method for computing the MPs from the CMB 
and LSS data. We also describe our simulations of CMB 
and LSS. In Appendix B we derive the analytical formula 
for the galaxy bispectrum. In Appendix C we compare the 
analytical MPs of CMB with non-Gaussian simulations in 
the Sachs- Wolfe limit. In Appendix D, we extend the cor- 
rections of primordial potential to n-th order, in order to 
examine more carefully validity of our perturbative expan- 
sion. 

Throughout this paper, we assume ACDM cosmology 
with the cosmological parameters at the maximum like- 
lihood peak from the WMAP flrst-year data only fit 
(Spergel et al. 2003). Specifically, we adopt = 0.049, 
^cdm = 0.271, Vtf, = 0.68, Ho = 68.2 km s"^ Mpc-\ 
T = 0.0987, and = 0.967. The amplitude of primor- 
dial fluctuations has been normalized by the first acoustic 
peak of the temperature power spectrum, 1{1+1)Ci/{2tt) = 
(74.7 jjKf at I = 220 (Page et al. 2003b). 

2. GENERAL PERTURBATIVE FORMULA FOR MINKOWSKI 

FUNCTIONALS 

Suppose that we have a d-dimensional fluctuating field, 
/, which has zero mean. Then, one may deflne the MPs 
for a given threshold, z/ = //(To, where (Tq = (./^)^^^ is the 
standard deviation of /. Matsubara (2003) has obtained 
the analytical formulae for the A;-th Minkowski Punction- 

als of weakly non-Gaussian fields in d-dimension, V^'^\i'), 
as (Eq. (133) of Matsubara 2003) 
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where Hn{i') are the Hermite polynomials, and cok = 
7r'^/^/r(fc/2 -f 1) gives u)o = 1, = 2, iJ2 = tt, and 
L03 = 47r/3. Here, 5^'^ are the "skewness parameters" 
defined by 
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Table 1 

Minkowski Functionals defined in d— dimensional space: d = 2 for CMB and d = 3 for LSS. 



observations 


d 


Vo 


Vi 


V2 


V3 


CMB 


2 


Area 


Total Circumference 


Euler Characteristic 




LSS 


3 


Volume 


Surface Area 


Total Mean Curvature 


Euler Characteristic 



qa)_ 3 (/^(vV)) 



2) 



3d ((V/j- (V/)(V-/)) 



(4) 



(5) 



2{d-l) o\ 

which characterize the skewness of fluctuating fields and 
their derivatives. The quantity Oi characterizes the vari- 
ance of fluctuating fields and their derivatives, and is given 
by 

a^^J^°"^e^P{k)W'{kR), (6) 
for d = 3, and 

'^]^^T.('^l + l)[l{l + i)YCiW^, (7) 

for d = 2. In both cases W represents a smoothing kernel, 

or a window function, which will be given by a product 
of the experimental beam transfer function, pixelization 
window function, and an extra Gaussian smoothing. The 
power spectra, P{k) for d = 3 and C; for d = 2, are defined 
as 

(/k/^') = (27r)3nfc)fe(k-k'), (8) 

(o-im «/*,„') = ClSu'Smm', (9) 

where Sd{^) is the Dirac delta function, and the Fourier 
and harmonic coefficients are given by 

/(x) = / 7— ^/ke^" •^ (10) 



(11) 



. (27r)3- 



Im 



for d = 3 and 2, respectively. Finally, the most relevant 
Hermite polynomials are given by 



^..M.;|e-/Wc(-^), 


(12) 


^^0 = 1, 


(13) 




(14) 




(15) 


H^{v)=v^ -3v, 


(16) 


Hi{v)=v'^ -Qv"^ +3, 


(17) 


Hsiu) =v^ - lOz/^ + 15i/. 


(18) 



3. APPLICATION l: COSMIC MICROWAVE BACKGROUND 

3.1. Analytical Formula for Minkowski Functionals of 

CMB 

For the cosmic microwave background, we have d = 2 
and / = AT/T. We define the angular bispcctrum as 

^lihh ^ ~ {^hmiai2m2'^l3m3) ■ (19) 



Then, by expanding skewness parameters into spherical 
harmonics, we obtain 



5(2) 



Lr.., 

3 ^ ll{h + l) + Hh + l) +13(13 + 1) 



E 



(21) 



8770-4 



E 

hrrii 



Il{ll + l) + l2{l2 + l)~l3{l3 + l)] 



xhih + 1) + (eye.)} BZ\i'"^^gZ7r' 

xWi.Wi^Wi,, (22) 

where (eye.) means the addition of terms with the same 
cyclic order of the subscripts as the previous term, Wi is 
a smoothing kernel in I space and Qi^l^t^"^^ is the Gaunt 
integral. 



jdn Yi,mA-^)Yi,mA-^)Yi^mA-^)- (23) 

Note that we have used the following properties of ^^(^(n): 
V^Ylm{n) = -l{l + l)Yim{n), (24) 

/ 



dh [VYi.mAi^)] . [VYi2m2{n)]Yi,mAi^) 

_ hih + l)+l2{l2+l)-l3{h + l) 

The summation over rrii can be done by writing 

-^/iWs =^hl2h °'l'2i3> (26) 

where h^i^i^^ is the reduced bispectrum that depends on 
specific non-Gaussian models (Komatsu & Spcrgcl 2001). 
Using the reduced bispectrum, we finally obtain the ana- 
lytical formula for MFs of the CMB: 

^^°^ = ^ E liM.Ki.hWi^\%Wi,, (27) 

2<h<l2<l3 



S^'^ = ^Z-2-2 E [h{h + l) + l2{l2 + l) + l3{h + l)] 

1 2<h<l2<l3 

xIiMsbhhhWi.Wi^Wi,, (28) 
3 



{Mh + l) + l2{l2 + l)-l3{h + l)] 



2<li<h<t3 



xhih + 1) + (eye.)} lf,i2i,biM3Wi,Wi2Wi,, (29) 



where 
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I1I2I3 



/(2Zi + l)(2Z2 + l)(2/3 + l) / h I2 I. 
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and we have used 

E/^jnim2m3\2 _ 7-2 /oix 

When /nl is a constant, the form of bi^i^i^ is given by 
(Komatsu & Spergel 2001) 

biM,=2 / r2dr[6/;(r)6f,(r)6i^^(r)+6/;(r)C''(r)6f3(r) 
Jo 

+''rW6f,(r)6f3(r)], (32) 

where 

6f (r) = - / k''dkP^{k)gTi{k)ji{kr), (33) 

6f^^(r) = - / k^dkfNLgTi{k)ji{kr), (34) 

and Pcj>{k) oc k"'~^ is the primordial power spectrum of 
The amphtude of ^^(fc) is fixed by the first peak ampli- 
tude of the temperature power spectrum, l{l+l)Ci/{27r) = 
(74.7 A«K)2 at I = 220 (Page et al. 2003a), and the tem- 
perature power spectrum is given by 

2 

Ci = - k^dkP^{k)g^i{k). (35) 

Jo 

We calculate the full radiation transfer function, griik), 
using the publicly-available CMBFAST code (Scljak & 
Zaldarriaga 1996). Note that our formalism is completely 
generic. One can easily generalize our results to non- 
Gaussian models with a momentum-dependent /nl by us- 
ing an appropriate form of bi-^i^i^ given in Liguori et al. 
(2006). Our results suggest that the MFs do not contain 
information beyond the bispectrum when non-Gaussianity 
is weak. The MFs of CMB basically measure the weighted 
sum of the CMB angular bispectrum. 

In Figure 1, we plot the skcwncss parameters, S^'^^ 
(Eqs. [27-29]), and S'^"^ multiplied by cto, for a pure sig- 
nal of CMB anisotropy (without noise) as a function of a 
Gaussian smoothing width, 9s, which determines a Gaus- 
sian smoothing kernel, Wi = cxp[— ^(Z + l)6'^/2]. The pcr- 
turbative expansion of MFs works only when S'^^Vq is 
much smaller than unity (see Eq. [2]). We find that the 
perturbative expansion is valid for /nl ^ 3300 from the 
results plotted in the right panel of Figure 1 that show 
|5('»)|ao <3x 10-4 /nl. 

This Figure also shows how MFs may be as powerful 
as the angular bispectrum in measuring /nl- Komatsu 
& Spergel (2001) have shown that sensitivity of the first 
skewness parameter, to /nl is much worse than that 
of the angular bispectrum, as acoustic oscillations in I 
space smear out non-Gaussian signals in the skewness, 
which is the weighted sum of the angular bispectrum over 
I. (The angular bispectrum is negative in the Sachs- Wolfe 
regime at low I, and oscillates about zero by changing its 
sign at higher I.) The MFs are sensitive to not only S^^\ 
but also S*^^-* and S^'^\ The weight of the sum over mul- 
tipoles differs among three skewness parameters: S'^^^ has 
the largest weight at high I, 5'-°^ has the largest weight at 
low /, and 5*^^^ is somewhere in between (Eqs. [27 29]). 
In particular, because S'^^-' picks up the highest multi- 
poles efficiently, S^^^ changes its sign depending on Og- 
5^^' is negative on very large angular scales. As Og de- 
creases (as the small scale information is included), 5^^^ 



increases, and eventually changes its sign to positive values 
near the scale of the first acoustic peak, 9s ~ 40 arcmin, 
where the bispectrum has the largest amplitude, while the 
other two skewness parameters do not change their signs. 
Therefore, S''^^^ keeps information about the acoustic os- 
cillations. This property is crucial for obtaining a better 
signal-to-noise ratio for primordial non-Gaussianity in the 
CMB. 

Figure 2 shows the predicted MFs of CMB temperature 
anisotropy, Vq, Vi and V2, as a function of /nl- The MFs 
for /nl = 100, 50 and are plotted in the solid lines, while 
the MFs for /nl = —50 and —100 are plotted in the dot- 
ted lines. The lower panels show the difference between the 
Gaussian and non-Gaussian MFs divided by the maximum 
amplitude of Gaussian MFs. Primordial non-Gaussianity 
with /nl = 100 changes Vq, Vi and V2 by 0.2%, 0.4%, 
and 1%, respectively, relative to the maximum amplitude 
of the corresponding Gaussian MFs. While Vq (area) has 
Uttle dependence on 9a, V2 (Euler characteristic) depends 
on 9s strongly, which is mainly due to the sign change of 
5(2) at 9s ~ 40 arcmin. 

In Appendix C we show that these analytical predictions 
agree with non-Gaussian simulations in the Sachs Wolfe 
limit very well. The comparison with the full simulations 
that include the full radiation transfer function will be re- 
ported elsewhere. 

3.2. Measuring /nl from Minkowski Functionals of CMB 

The MFs at different threshold values, u, are correlated, 
and different MFs arc also correlated, i.e., 

EitU^, ^ {vt'\'^)V^'\'^')) ^ ~ v'). (36) 

where (5 denotes the Kronecker delta. Therefore, it is im- 
portant to use the full covariance matrix, E, in the data 
analysis. We obtain the full covariance matrix from Monte 
Carlo simulations of Gaussian temperature anisotropy. 
Because non-Gaussianity is weak, the covariance matrix 
estimated from Gaussian simulations is a good approxi- 
mation of the exact one. In Appendix A.l we describe 
our methods for computing MFs from the pixelized CMB 
maps. In Appendix A. 2 and A. 3 we describe our methods 
for simulating sky maps of CMB temperature anisotropy 
with noise and instrumental characteristics of the WMAP 
and Planck experiment respectively. 

We use the Fisher information matrix formalism to esti- 
mate the projected errors on /nl from given measurement 
errors on MFs. The Fisher matrix, Fij , is written in terms 
of the inverse of the covariance matrix, S"^, as 

kk' 

(37) 

where pi is the i-th parameter. For CMB, we consider 
only one parameter, pi = fML, whereas for LSS we also 
include galaxy biasing parameters (see § 4). The projected 
l-cT error on /jvl is given by the square root of {F~^)ii. 
While equation (37) may be evaluated for a given smooth- 
ing scale, 9s, one will eventually need to combine all combi- 
nations of MFs at different 9s to obtain the best constraint 
on /nl from given data. The MFs at different 6s are also 
correlated. We therefore calculate the full covariance ma- 
trix of MFs that consists of N^, x Nmf x -^s elements, where 
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Fig. 1. — (Left) Skewness parameters, S^°'\ for a = (solid), 1 (dotted), and 2 (dashed). (Right) Skewness parameters multiplied by 
variance, S'^'cro, for a = (solid), 1 (dotted), and 2 (dashed). Both have been divided by /nl, and are plotted as a function of a Gaussian 
smoothing width, 6s. Note that S'^' changes its sign at 9s ~ 40 arcmin. 
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Fig. 2. — Analytical predictions for the Minkowski Functionals of CMB temperature anisotropy with primordial non-Gaussianity character- 
ized by /nl = —100, —50 (dotted), 0, 50, and 100 (solid). Each MF, Vj. (k = 0, 1, and 2), is plotted in the top panels. The other panels show 
the difference between non-Gaussian and Gaussian MFs, , divided by the maximum amplitude of V^*^. From the top to bottom, 9s = 10, 
20, 40, 70, and 100 arcmin are shown. 



= 25 is the number of bins for v per each MF in the 



range of i' from —3 to 3, A^mf = 3 and 4 is the number of 
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Table 2 

Projected l-tr errors on /nl from each of the Minkowski Functionals of CMB temperature anisotropy (Vq: 

AREA, Vi: total CIRCUMFERENCE, V2 : EULER CHARACTERISTIC) AS WELL AS FROM THE COMBINED ANALYSIS OF ALL THE 

Minkowski Functionals, for various combinations of smoothing scales, 9s. We use noise and beam properties of 
WMAP 1-YEAR and 8-year observations with the KpO mask. The last block of Table shows the 1-a errors from 
the "ultimate WMAP", which uses all sky with zero noise. (The beam smearing is still included.) 



WMAP 1-year WMAP 8-year no noise and no sky cut 



^ a L'^-'j-^ ^^^^^^^^ 


All 


{Vo,Vi,V2) 


All 




All 




100 


271 


(671,297,494) 


269 


(668,296,494) 


135 


(461,146,182) 


70 


153 


(456, 168, 272) 


153 


(456,168,271) 


94 


(363, 100, 137) 


40 


86 


(283,96,148) 


84 


(289,95,144) 


59 


(236,67, 94) 


20 


54 


(176,67,82) 


51 


(183,63,79) 


41 


(160,52,62) 


10 


73 


(139,92,91) 


41 


(134, 54,61) 


35 


(120,47,51) 


5 






75 


(125,94,87) 


30 


(106,42,43) 


10,20 & 40 


46 


(92,62,67) 


36 


(82,51,56) 


30 


(49,43,45) 


5, 10 & 20 






36 


(71,52,53) 


22 


(46, 32, 37) 



Table 3 

Same as Table 2 but for Planck data. 



or ■ 1 Planck no noise and no sky cut 

(^^arcmmj (Up, t^i, K2) AH (t/p, Ki, K2) 

100 27l (665,299,502) l36 (462, 147, 184) 

70 157 (463,171,278) 95 (354,102,141) 

40 77 (281,86,128) 54 (235,61,84) 

20 46 (153,59,69) 37 (140,49,55) 

10 32 (115,45,45) 28 (102,40,39) 

5 24 (97,31,36) 21 (87,27,30) 

5,10fc20 19 (46,26,31) 15 (38,22,26) 



MFs for = 2 and 3, respectively (i.e., Nmf = 3 for CMB 
and 4 for LSS), and Ns is the number of smoothing scales 
used in the analysis. The Fisher matrix may be written as 

-2.^(^ ^^^^ 

aa' 

where a is a single index denoting fc, and 0^- 

Let us comment on the effect of noise on MFs. The 
instrumental noise increases cr? by adding extra power at 
small scales. On the other hand, the signal part of the 
angular bispectrum is unaffected by noise because noise 
is Gaussian. As a result, the instrumental noise always 
reduces the skewness parameters, as the skcwness param- 
eters contain aj in their denominator. Therefore, the 
MFs would approach Gaussian predictions in the noise- 
dominated limit, as expected. We compute the increase 
in cTp and cti due to noise from Monte Carlo simulations, 
and then rescale 5^"' and dVa / 9/nl for a given smoothing 
scale, Og. The window function, Wi, includes the beam 
smearing effect, pixel window function, and a Gaussian 
smoothing. 

Using the method described above and in Appendix A, 
we estimate the projected 1-a error on /jvl expected from 

WMAP's 1-year and 8-ycar observations. We also con- 
sider an ideal WMAP experiment without noise or sky 

^ Note that 0s is a Gaussian width, which is 1/2.35 times the full w 



cut (but the beam smearing is still included). We consider 
six different smoothing scales, 9s — 5,10,20,40,70, and 
100 arcmin. The results from various combinations of Og 
are summarized in Table 2 for WMAP observations. 

For the WMAP 1-year data, the error on /nl is the 
smallest for 9s = 20 arc-minutes. At smaller angular 
scales, say 9s = 10 arc-minutes, the noise dominates more 
and thus the error on /nl increases at 9s ^ 20 arc-minutes. 
For the WMAP 8-year data, on the other hand, a better 
signal-to-noise ratio at smaller angular scales enables us to 
constrain MFs at 9s = 10 arc-minutes. As the beam size 
of WMAP in W band is about 9s — 10 arc- minutes'^, one 
cannot constrain MFs at the angular scales smaller than 
this. When all the smoothing scales are combined, the pro- 
jected 1-a error on /nl reaches ~ 40 for the WMAP data, 
which is in a rough agreement with the result reported in 
Komatsu et al. (2003) and Spergel et al. (2006). The best 
constraint that can be obtained from the WMAP data, in 
the limit of zero noise and full sky coverage, is /nl ^ 22. 

We also estimate the Planck constraint on /nl listed at 
Table 3. As Planck's beam and noise are ~ 4 and 10 times 
as small as WMAP's, respectively, one can constrain MFs 
even at = 5 arc-minutes. Planck should be sensitive to 
|/nl| ~ 20. 

at half maximum. 
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4. APPLICATION II: LARGE-SCALE STRUCTURE 

4.1. Analytical Formula for Minkowski Functionals of 

LSS 

For the large-scale structure, we have d = 3 and / = 
(5g(x, z), where dg is the density contrast of galaxies. Sta- 
tistical isotropy of the universe gives the following form of 
the bispectrum: 

(4(0)(ki)4(z)(k2)4(2)(k3)) = {2TrfdD0^i + k2 + ka) 

xBg{ki,k2,k3,z). (39) 
We obtain the skewness parameters by integrating 
Bg{k\,k2,kz,z) over fci, k2, and ji = (ki ■ \i.2) / {kik^) with 
appropriate weights as 

1 



dki 



pOO pi 

L */. 



5(1) 

g 



{z) = 



xfc^fc|i?g(fci, A:2, fci2, 
X W{kiR)W(k2R)W{ki2R) 
1 



u2 7^2 



/ dki / rffca / 
Jo Jo J-i 



xkik2{ki + ^2 + nkik2)Bg{ki,k2, ki2, z) 
xW{kiR)W{k2R)W{ki2R), 
3 



/■OO /-OO /■! 

/ dki / (ifc2 / c?/i 
Jo Jo J -I 



mn^al^{z) 
xkfk^il - ii^)Bg{ki,k2, ki2, z) 
xW{kiR)W{k2R)W{ki2R), 



(40) 
dfj. 

(41) 
(42) 



where fci2 = |ki + k2| = (fc 



^2 



U/2 



Unlike for the CMB, where we needed to consider only 
the effect of primordial non-Gaussianity, there are three 
sources of non-Gaussianity in Bg(ki,k2, k3,z): primordial 
non-Gaussianity, non-linearity in the gravitational evolu- 
tion, and non-linearity in the galaxy bias. In Appendix B 
we show that Bg is given by 

Bg{ki,k2,ks,z) = bl{z) [BpHm{ki,k2,k3,z) 

+Bgra,-v{ki, k2, ks, z)] + bl{z)b2{z) 
X [Pm(fci, 0)P^(fc2, z) + (eye.)] , (43) 
where Pm{k, z) is the linear matter power spectrum, bi{z) 
and b2{z) are the linear and non-linear galaxy bias param- 
eters, respectively (sec Eq. [Bl] for the precise definition), 
and -Bprim and .Bgrav represent the contributions from pri- 
mordial non-Gaussianity and non-linearity in gravitational 
clustering, respectively: 



-Bpriin(A;i,fc2,fc3,^) = 



_ 2/nl 



Pm{ki,z)Prn{k2,z)M{k3) 

M{ki)M{k2) 



D{z) 

+(cyc.)] , (44) 

Bgr:scv{k\,k2,k3,z) = 2 [F2 (ki , k2)Pm (^1 , ^;)Pin(fc2 , 2) 

+(cyc.)] , (45) 

where D{z) is the growth rate of linear density fluctua- 
tions normalized such that D{z) — > 1/(1 -I- z) during the 
matter era, and M{k) and P2(ki,k2) are given by equa- 
tion (BIO) and (B5), respectively. These equations suggest 
that /nl and the galaxy bias parameters must be deter- 
mined simultaneously from the LSS data. Moreover, even 
if the galaxy bias is perfectly linear, 62 = 0, the primor- 
dial signal might be swamped by non-Gaussianity due to 
non- linear gravitational clustering, Pgrav 



In order to investigate how important the effect of -Bgrav 
is, let us define the skewness parameters that are con- 
tributed solely by -Bprim or .Bgrav Substituting i^grav and 
(Jm,j{z) for -Bg and Ugj, respectively, in equation (40), 

(41), and (42), we obtain Sgrlv (a = 0, 1, and 2) given 
by (Matsubara 2003) 



'^erav('^) ^ 



287rV4,o(z) 
3 



[5/220 (-2) + 7/131 (-2) + 2/222(^)1 , 

(46) 

[10/240 (;j) + 12/331 (Z) 



9(1) (z\ = 

''^■■-^'^^ 56^V2,,,(z)<^,(z) 

+7/151 {z) + 1 1/242 {z) + 2/333 [z)], (47) 



'^srav('^) ~ 



9 



^ 567rV4,i(^) 

-7/353(^) - 2/444(^)1 , 



[5/440 (-Z) + 7/351(2) - 3/442 (-2) 

(48) 



where 



-(z) = / dk 
Jo 



if dk2 [ dfiW{kiR)W(k2R) 
Jo J-i 



:W{ki2R)kTk^ti^P^{ki,z)P^{k2,z). (49) 



Note that 



61 (T 



Similarly, we also calculate 



the primordial skewness parameters, S^l^, by substitut- 
ing -Bprim and (Jni.j(^) foi' -Bg and o-gj{z), respectively, in 
equation (40), (41), and (42). These skewness parameters 
are related to the skewness parameters of the total galaxy 
bispectrum as 



■ 



^(a) 
prim 



I c(o) 
~r ^grav 



bl 



+ 



362 
6?- 



(50) 



In Figure 3 we compare S^l^{z) and S'grav(^) at z = 
0, 2, and 5, as a function of a smoothing length, /?, 
which is in units of h^^ Mpc. (The smoothing ker- 
nel W{kR) is set to be a Gaussian filter, W{kR) = 
exp(— fc^/i^/2).) Non-Gaussianity from non-linear gravi- 
tational clustering always gives positively skewed density 
fluctuations, S'grlv > 0. Primordial non-Gaussianity with 
a positive /nl also yields positively skewed density fluctu- 
ations; however, primordial non-Gaussianity with a neg- 
ative /nl yields negatively skewed density fluctuations, 
which may be distinguished from S'|rav more easily. As the 
smoothing scale, R, increases (i.e., density fluctuations be- 
come more linear), non-Gaussianity from non-linear clus- 
tering, S'gri,vO'm,o, becomes weaker, while primordial non- 
Gaussianity, S^l^a-m.o, remains nearly the same. At z = 
0, the primordial contribution exceeds non-linear gravity 
only at very large scales, R > 200 Mpc for /nl = 100, 
and R > 800 Mpc for /nl = 10. As higher redshift, on 
the other hand, non-linearity is much weaker and therefore 
the primordial contribution dominates at relatively smaller 
spatial scales, R > 120 Mpc and 80 Mpc at z = 2 
and 5, respectively, for /nl = 100. Unlike for the CMB, all 
the skewness parameters of galaxies exhibit similar depen- 
dence on the smoothing scales. The perturbation formula 
is valid when the amplitude of the second order correction 
of MFs is small, S^'^^ao « 1, that is, /nl < 5000. 
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R[h-'Mpc] R[h-'Mpc] 



Fig. 3. — (Left) Skewness parameters, S^°'"'{z), for a = (solid), 1 (dotted), and 2 (dashed), at z = 0, 2, 5. (Right) Skewness parameters 
multiplied variance, S'(°)(z)(To, for a = (solid), 1 (dotted), and 2 (dashed), at z = 0, 2, 5. In both fi Hures, the thick lines show the skewness 
parameters from non-linear gravitational clustering, Sgrlv (Eq. [46-48]), while the thin lines show the primordial skewness parameters, S^^)^^^^, 
with /nl = 100 (and 10 in the right panel). The symbols show the skewness parameters measured from numerical simulations of non-Gaussian 
matter density fluctuations. Note that Sgrlv and s'^^^^ao are independent of z. 



4.2. Measuring /nl from Minkowski Functionals of LSS 

Figure 4 shows the perturbation predictions for the MFs 
from primordial non-Gaussianity with /nl — 100, 50, 0, 
—50, and —100. For comparison, we also show the MFs 
computed from numerical simulations with /nl = 100. 
In Appendix A. 4 we describe our methods for computing 
MFs from the LSS data. In Appendix A. 5 we describe 
our methods for simulating the LSS data with primordial 
non-Gaussianity (but without any effects from non-linear 
gravitational clustering or galaxy bias). The error-bars are 
esti mated from variance among 2000 realizations divided 
by V2000. The left panels show the MFs, while the right 
panels show the difference between the non-Gaussian and 
Gaussian MFs, divided by the maximum amplitude of each 
MF. We find that the analytical perturbation predictions 
agree with the numerical simulations very well. 

Let us comment on some subtlety that exists in the 
comparison between the perturbation predictions and nu- 
merical simulations. The MFs measured from numerical 
simulations often deviate from the analytical predictions, 
even for Gaussian fluctuations, due to subtle pixelization 
effects. The MFs from our Gaussian simulations deviate 
from the analytical predictions at the level of 10~^ when 
normalized to the maximum amplitude of each MF. It is 
important to remove this bias, as the magnitude of this ef- 
fect is comparable to or larger than the effect of primordial 
non-Gaussianity with /nl = 100 (lO^'^ for Vq and 10^^ 
for V3). Therefore, it is often necessary to re-calibrate the 
Gaussian predictions for the pixelization effects by run- 
ning a large number of Gaussian realizations. However, 
we have found that the difference between the Gaussian 
and non-Gaussian MFs measured from simulations does 
agree with the perturbation predictions without any cor- 



rections. Therefore, one may use the following procedure 
for calculating the correct non-Gaussian MFs: 

(1) Use the analytical formulae (Eq. [2]) to calculate the 
difference between the Gaussian and non-Gaussian 
MFs, 

AV^''\'y, /nl) ^ V^'\,y, /nl) - V^'\,,, /nl = 0). 

(51) 

(2) Run Gaussian simulations. Estimate the average 

MFs from these Gaussian simulations, V'^f'*'' {v, /nl = 
0) . This would be slightly different from the analyti- 
cal formula, V^'^'^ (v, /nl ~ 0), due to the pixelization 
and boundary effects. Note that the same simula- 
tions may be used to obtain the covariance matrix 
of MFs. 

(3) Calculate the final non-Gaussian predictions as 

Vt''\y, /nl) = vt'Hu, /nl = 0) + ^Vt''\y, /nl). 

(52) 

We estimate the projected errors on /nl from LSS using 
the Fisher information matrix in the same way as CMB 
(see § 3.2). We compute the covariance matrix of MFs 
from 2000 realizations of Gaussian density fluctuations 
in a 1 Gpc'^ cubic box. For simplicity, we focus on the 
large scales and ignore non-Gaussianity from non-linear 
gravitational evolution. (The minimum smoothing scale is 
R — QO Mpc). We also ignore shot noise. 

Table 4 shows the projected 1-cr errors on /nl from 
a galaxy survey covering 1 Gpc'^ volume. This vol- 
ume would correspond to e.g., a galaxy survey covering 
300 deg^ on the sky at 5.5 < z < 6.5. Note that these 
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Fig. 4. — Analytical predictions for the Minkowski Punctionals of LSS with primordial non-Gaussianity characterized by /nl = 100, 50, 
(solid), —50, and —100 (dotted), for a Gaussian smoothing length of R = 100 Mpc. The left panels show the MFs Vj; (fc = 0, 1, 2, and 
3 from top to bottom), while the right panels show the difference between non-Gaussian and Gaussian MFs, V^, divided by the maximum 
amplitude of . The symbols show the average and error of the MFs calculated from 2000 realizations of simulated primordial non-Gaussian 
density fluctuations with /nl = 100. 
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constraints are independent of z when non-linear gravi- 
tational evolution is ignored, as S^l^ao is independent 
of z. We have assumed a linear galaxy bias (62 ~ 0) in 
the second column, while we have marginalized 62/^1 in 
the third column. A more realistic prediction would lie 
between these two cases, as we can use the power spec- 
trum and bispectrum to put some constraints on 62- We 
find bettor constraints from smaller smoothing scales for a 
given survey volume, simply because we have more modes 
on smaller scales. The limits on /nl from a 1 Gpc"* 
survey are not very promising, |/nl| ^ 270 at the 68% 
confidence level; thus, one would need the survey volume 
as large as 25 /i^"^ Gpc"^ to make the LSS constraints com- 
parable to the WMAP constraints (using the MFs only). 
This could be done by a survey covering ^ 2000 deg^ at 
3.5 < z < 6.5. One would obviously need more volume to 
make it comparable to the Planck data. 

5. SUMMARY AND CONCLUSIONS 

We have derived analytical formulae of the MFs for 
CMB and LSS using a perturbation approach. The analyt- 
ical formula is useful for studying the behavior of MFs and 
estimating the observational constraints on /nl without 
relying on non-Gaussian numerical simulations. The per- 
turbation approach works when the skewness parameters 
multiplied by variance, S'^ctq, is much smaller than unity, 
i.e., |/nl| < 3300 for CMB and |/nl| < 5000 for LSS, 
both of which are satisfied by the current observational 
constraints from WMAP (Komatsu ct al. 2003; Crcminelli 
et al. 2006; Spergel et al. 2006). We have shown that the 
perturbation predictions agree with non-Gaussian numer- 
ical realizations very well. 

We have used the Fisher matrix analysis to estimate the 
projected constraints on /nl expected from the observa- 
tions of CMB and LSS. We have found that the projected 
1-(T error on /nl from the WMAP should reach 50, which 
is consistent with the MF analysis given in Komatsu et 
al. (2003); Spergel et al. (2006), and is comparable to the 
current constraints from the bispectrum analysis given in 
Komatsu et al. (2003); Creminelh et al. (2006); Spergel et 
al. (2006). The MFs from the WMAP 8-year and Planck 
observations should be sensitive to |/nl| ~ 40 and 20, re- 
spectively, at the 68% confidence level. 

As the MFs arc solely determined by the weighted 
sum of the bispectrum for |/nl| ^ 3300 for CMB and 
|/nl| ^ 5000 for LSS, the MFs do not contain informa- 
tion more than the bispectrum. However, this docs not 
imply that the MFs are useless for measuring primordial 
non-Gaussianity by any means. The important distinction 



between the MFs and bispectrum is that the MFs are in- 
trinsically defined in real space, while the bispectrum is 
defined in Fourier space. The systematics in the data are 
most easily dealt with in real space, and thus the MFs 
should be quite useful in this regard. Therefore, in the 
presence of real- world issues such as inhomogeneous noise, 
foreground, masks, etc., these two approaches should be 
used to check for consistency of the results. 

In this paper we have calculated the MFs from primor- 
dial non-Gaussianity with a scale-independent /nl- It is 
easy to extend our calculations to a scale-dependent /nl- 
All one needs to do is to calculate the form of the bispec- 
trum with a scale-dependent /nl (e.g., Babich et al. 2004; 
Liguori et al. 2006), and use it to obtain the skewness 
parameters, S'-"^ (Eqs. [20-22] for CMB and Eqs. [40- 
42] for LSS). The MFs are then given by equation (2) 
in terms of the skewness parameters. Also, we have 
not included non-Gaussianity from secondary anisotropy 
such as the Sunyaev-Zel'dovich effect, Rees-Sciama effect, 
patchy reionization, weak lensing effect, extragalactic ra- 
dio sources, etc. It is again straightforward to calculate 
the MFs from these sources using our formalism, as long 
as the form of the bispectrum is known (e.g., Spergel & 
Goldberg 1999; Goldberg & Spergel 1999; Cooray & Hu 
2000; Komatsu & Spergel 2001; Verde & Spergel 2002). 

The constraints on primordial non-Gaussianity from the 
MFs of LSS in a galaxy survey covering 1 h^^ Gpc"^ vol- 
ume are about 5 times weaker than those from the MFs of 
CMB in the WMAP data. One would therefore need the 
survey volume as large as 25 Gpc^ to make the LSS 
constraints comparable to the WMAP constraints (using 
the MFs only). This could be done by a survey covering 
^ 2000 deg^ at 3.5 < z < 6.5. One would obviously need 
more volume to make it comparable to the Planck data. 
The MFs from LSS are less sensitive to primordial non- 
Gaussianity because non-Gaussianity from the non-linear 
evolution of gravitational clustering exceeds the primor- 
dial contribution at i? < 200, 120, and 80 h-^ Mpc at 
^; = 0, 2, and 5, respectively, which severely Hmits the 
amount of LSS data available to constrain primordial non- 
Gaussianity. 
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Constraints on /nl from the MFs of LSS in a galaxy survey covering 1 Gpc^ volume. All of the four MFs 

HAVE BEEN COMBINED. Tl-IE FIRST COLUMN SHOWS THE SMOOTHING SCALES THAT HAVE BEEN USED. THE SECOND COLUMN 
SHOWS THE CONSTRAINTS ASSUMING A LINEAR BIAS (62 = 0), WHILE THE THIRD COLUMN SHOWS THE CONSTRAINTS WITH 62/61 
MARGINALIZED. NOTE THAT THESE CONSTRAINTS ARE INDEPENDENT OF z WHEN NON-LINEAR GRAVITATIONAL EVOLUTION IS 

IGNORED, AS S^ll^ao IS INDEPENDENT OF Z. 
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APPENDIX 

A. MEASURING MINKOWSKI FUNCTIONALS FROM CMB AND LSS 

In this Appendix we describe our methods for computing the MFs from the CMB and LSS data. We also describe our 
simulations of CMB and LSS. 

A.l. Computational Method: CMB 

We estimate the MFs from pixelized CMB sky maps by integrating a combination of first and second angular derivatives 
of temperature anisotropy, Ik, over the sky (Schmalzing & Gorski 1998), 



V^f(-) = ^^|^#^, (Al) 



where flj is the unit vector pointing toward a given position on the sky. We set the weight of i—th pixel, Wi, to be 1 
when the pixel at fli is outside of the survey mask, and otherwise. We calculate Tk at flj from covariant derivatives of 
temperature anisotropy divided by its standard deviation, u{fti) = {AT/T)/ao, 

To = e{u-v), (A2) 

I, = If{u-,.)^u^, + u%, (A3) 

I, = -Fiu-.) . (A4) 

The function F{u — v) has a value of 1/Ai/ (Ai/ is the binning width of v) when u is within \v — AZ//2, v + Av/2], and 
otherwise. Covariant derivatives are related to the partial derivatives as 

u-0 = ufi, (A5) 

U:4> = -^u,^, (A6) 
sm 

u-00 = Ufi0, (A7) 
1 cos^ 

u-0<j> = -T—p:U^0^ r-o-^u^^,, (A8) 

sin^ ' sin^ 

1 cos9 
sm smti 
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Table A5 

Sky fraction used for the calculation of MFs 



9s sky fraction [%J 



100 


13 


70 


26 


40 


41 


20 


62 


10 


73 


5 


75 



The first derivative of the temperature field u^i(i = 0, or (f)) are calculated in Fourier space as; 



(AlO) 



yim,e 



I 



tan^ 



(m = /) 



The second derivatives are also computed as 



m,ij 



Ira 



Yim,ee 



Yim,e<t, = imYimfi 



Ylr. 



tan^ 



(All) 

(A12) 
(A13) 

(A14) 

(A15) 

(A16) 
(A17) 

We calculate the MFs of temperature anisotropy for u from —3 to 3 with Aj, = 0.24, which yields 25 bins per each MF. 

A. 2. Simulating WMAP Data 

In order to quantify uncertainties in the estimated MFs from cosmic variance and various effects, we simulate Gaussian 
temperature anisotropy maps with noise characteristics of the WMAP data. We generate 1000 Gaussian realizations of 
CMB temperature anisotropy using the HEALPix package (Gorski et al. 2005). We use A^sido = 128, 256 and 512 for 
9s > 40, 40 > 6's > 10, 10 > 9s respectively. The number of pixels is given by ripix = 12N^^^^. From each sky reahzation 
we construct eight simulated maps of WMAP differential assemblies (DAs), Ql, Q2, VI, V2, Wl, W2, W3, and W4, by 
convolving the sky map with the beam transfer function in each DA (Page ct al. 2003a), and adding independent Gaussian 
noise realizations following the noise pattern in each DA. Each pixel is given noise variance of ctq ^^isi^/Nobsi^i): where 
-^obs(^i) is the number of observations per pixel, and ctq, noise given in Bennett et al. (2003a). We then co-add the 
eight maps by weighting each map by noise where TVo^s is the full-sky average of A^obs(^i)- Finally, we mask the 

co-added map by the KpO Galaxy mask (including point-source mask) provided by Bennett et al. (2003b). This mask 
leaves 76.8% of the sky available for the subsequent data analysis. In addition to WMAP one-year data, we also simulate 
the future WMAP eight-year data, by simply multiplying A^obs by a factor of 8. 

Before we estimate the MFs from each simulated map, we smooth it using a Gaussian filter with a smoothing scale of 



Wl = exp 



(A18) 



To remove the effect of survey mask, we calculate the Minkowski Functionals by limiting the pixels where the five 
measurements, cr,ai, and S'''^\k = 0,1, &2), have nearly equal values between the field with survey mask and that 
without survey mask. We use the pixels where the difference is within 5% of each standard deviation for 9s = 20 and 40 
arcmin and 10% for 9s = 70 and 100 arcmin. For 9s = 5 and 10 arcmin, we use the pixels as long as they arc away from 
the boundary of the mask by more than 29s- Table A5 lists the sky fraction used in the analysis for each smoothing scale. 
The mean density over the pixels which are used for the calculation of MFs is not completely zero and thus we subtract 
it from the density field to satisfy zero mean in each realization. 
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Table A6 

Approximate instrument specifications of the Planck satellite (http://www.rssd.esa.int/Planck). The central 

FREQUENCY, FULL WIDTH AT HALF MAXIMUM (FWHM) OF BEAM SIZE, 6Ifwhm, AND NOISE PER PIXEL (AT/T PER ^|whm) FOR 
14 MONTHS OF OBSERVATIONS ARE SHOWN PER EACH FREQUENCY CHANNEL. 



Instrument 




LFI 










HFI 






Center Frequency [GHzJ 


30 


44 


70 


100 


143 


217 


353 


545 


857 


Beam size ^fwhm [arcmin] 


33 


24 


14 


9.5 


7.1 


5 


5 


5 


5 


Pixel Noise AT/T [10*^] 


2 


2.7 


4.7 


2.5 


2.2 


4.8 


14.7 


147 


6700 



A. 3. Simulating Planck Data 

For the simulations of the Planck data, we follow the same procedure as the WMAP simulations. Each realization 

is a coaddcd map of 9 bands with the inverse weight of the noise variance listed in Table A6. Wc approximate the 
beam transfer function as a Gaussian function, exp{—9^^^^l{l + l)/2), where ^beam = ^fwhm/ v8hi(2). We add the 
homogeneous noise distribution with the noise variance per pixel cr^^j^^ given by 

noise i \ / I'^FWHmJ 

where i denotes each band, and ripix = ISiV^jj^, is the number of pixels in simulated maps. Wc use the KpO mask to define 
the survey area for the Planck simulations, in exactly the same manner as for the WMAP simulations. 

A. 4. Computational Method: LSS 

We use two complementary routines to compute the MFs of a density field on the grids. The first approach is often 
called Koenderink invariants (Koenderink 1984) in which the surface integrals of the curvature are transformed into the 
volume integral of invariants formed from the first and second derivatives of the density fluctuations. The second method, 
which is called Crofton's formula (Crofton 1868; Schmalzing & Buchert 1997; Koenderink 1984), is based on the integral 
geometry and the calculation reduce to simply counting the elementary cells (e.g., cubes, squares, lines, and points for 
the cubic meshes). The outline of these methods are summarized in Schmalzing & Buchert (1997) and the observational 
application to SDSS galaxy samples are performed by Hikage et al. (2003). 

A. 5. Simulating LSS Data with Primordial Non-Gaussianity 

We calculate the MFs of density field in a cubic box with a length of 1 Gpc. but ignore the observational effects 
such as survey geometry, for simplicity. We also ignore non-linear gravitational clustering or galaxy bias in order to isolate 
the effect from primordial non-Gaussianity. (The purpose of this simulation is to check the accuracy of our perturbation 
predictions for the form of MFs from primordial non-Gaussianity.) To simulate the LSS data with primordial non- 
Gaussianity, we first generate a Gaussian potential field in a cubic box with a length of 1 Gpc, assuming that the 
power spectrum of potential is P^pik) cx k""^^ where rig = 0.967. We inversely Fourier-transform it into real space to obtain 
0(x). (The number of grids is 128^.) We then construct a non-Gaussian potential field, <f>(x), using equation (eq.[l]) for 
a given /nl- We finally convert it to the matter density field by multiplying $ by M{k) in Fourier space (see Eq. [BIO]). 
We have generated 2,000 realizations of the non-Gaussian density field. 

B. DERIVATION OF GALAXY BISPECTRUM 

In this Appendix we derive the perturbative formula for the galaxy bispectrum including primordial non-Gaussianity, 
non-linearity in gravitational clustering, and non-linearity in galaxy biasing, in the weakly non-linear regime (Verde et al. 
2000; Scoccimarro et al. 2004). In the weakly non-linear regime, it would be reasonable to assume that the galaxy biasing 
is local and deterministic. We then expand the galaxy density contrast, Sg, perturbatively in terms of the underlying 
matter density contrast, dm, as (Fry & Gaztanaga 1993) 

S,iz) = bo{z) + h{z)6Uz) + ^-^5l{z) + 0(e). (Bl) 

where b^iz) is determined such that ((5g(z)) = 0. Here, hi{z) and 62 (^) are the time-dependent galaxy bias parameters. 
The power spectrum and bispectrum of the galaxy distribution, Pg and Bg, respectively, are then given by those of the 
underlying matter distribution as 

P^{k,z)^hl{z)PUk,z), (B2) 

Sg(fci, fc2, fcs, z) = b\{z)Bm{ki,k2, fcg, z) + hl{z)h2{z)[Pm{ki,z)Pm{k2, z) + (cyc.)], (B3) 
respectively. If the underlying mass distribution obeyed Gaussian statistics, its bispectrum would vanish exactly, 5^ = 0; 
however, the non-linear evolution of density fluctuations due to gravitational instability makes (Jm slightly non-Gaussian 
in the weakly non-Gaussian regime, yielding non-zero bispectrum. 
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The second-order correction to the density fluctuations from non-Hnear gravitational clustering gives the following 
equation, 

L,k(2:) = ^L,k(2:) + j d^qF2{(l, k - q)5L,q(-2)^L,k-q(-2), (B4) 

where Si,^k{z) is the linear (but non-Gaussian) density fluctuations, and 



F.(k„k,) = ^ + ^(^ + ^)+;-i^, (B5) 



5 ki ka (ki kA 2 (ki ■ 
7 2fcifc2 \k2 kj 7 klkl 
is the timc-indcpcndent kernel describing mode-coupling due to non-linear clustering of matter density fluctuations in 
the weakly non-linear regime. Equation (B5) is exact only in an Einstein-de Sitter universe, but the corrections in other 
cosmological models are small (e.g., Bernardeau 1994). The power spectrum and bispectrum of the underlying mass 
density distribution, and i?m, are thus given in terms of the linear and non- linear contributions: 

P^{k,z) = Pi^{k,z), (B6) 
Bra{kuk2, ks, z) = Bpri,„(fci, A:2, k^,z) + 2 [F2{\<iiM)Pi.{ku z)Pi^{k2, z) + (eye.)] . (B7) 
Note that we have ignored the non-linear contributions in the power spectrum. That is to say, the power spectrum is still 
described by linear perturbation theory. 

The remaining task is to relate 5i,{z) to Bardeen's curvature perturbations during the matter era, One may use 
Poisson's equation for doing this: 

fc'^kT(fc) = At:Gp^{z)^0^ = ^a„i/o'^m,k(^)(l + z), (B8) 

where T{k) is the linear transfer function that describes the evolution of density fluctuations during the radiation era and 
the interactions between photons and baryons (Eisenstein & Hu 1999). Note that $ is independent of time during the 
matter era. At very early times, say, z = z^ 10^, the non-linear evolution may be safely ignored at the scales of interest, 
and one obtains 

~ M(fc)dk 

0L,k(2;*) = -p- , (B9) 

1 + -2* 

where 

Therefore, using the quadratic non-Gaussian model given in equation (1), one obtains the linear power spectrum at ^ = 2;,, 

[l + z^ 

and the primordial bispectrum aA, z = z^, 

M{k i)M{k2)M {k3) , 
{l + z.f 

'P^{kuZ,)Pl^{k2,Z,)M{k^) 



= n , 1 L P't'ik), (BU) 



Bprin.{ki,k2,ks,z.) = 2/nl ^ , ! [" [P<t>iki) P4,{k2) + (eye.)] 



2/nl(1 + z*) 



(eye.) 



(B12) 



M{ki)M{k2) 

where P,^(A;) oc fc"""^ is the power spectrum of (j) and we have ignored the higher-order terms. We then use the linear 
growth rate of density fluctuations, (5l oc .0(2;), to evolve the linear bispectrum forward in time: 

'Pi.{k^,z)Pi^{k2,z)M{k^) 



-Dprim(Kl, K2, k3,Z) = 2/nL- 



M{ki)M{k2) 



(eye.) 



(B13) 



D{z) 

One may simplify this expression by normalizing the growth rate such that 

D{z,) = (B14) 

Note that this condition gives the normalization of D{z) that is actually independent of the choice of z», when is taken 
to be during the matter era. 

By putting all the terms together, we finally obtain the following form of Bg{ki, fc2, /ca, z): 

'Pm{ki,z)Pn,{k2,z)M{k3) 



Sg(fcl,fc2,fc3,^) = 2/NL^ 



+ (eye.) 



M{ki)M{k2) 
+2bl{z) [F2{ki,-k2)Pm{ki,z)P^{k2,z) + (eye.)] 

+bl{z)b2{z)[Pn,{ki,z)P^{k2,z) + (eye.)]. (B15) 

We use this formula to calculate the skewness parameters that are used for the MFs of the galaxy distribution. 

While the perturbative formula for the MFs derived by Matsubara (2003) (See § 2) works well for i>, it is not immediately 
clear if it works for Sg because of the /c-dependent coeflacient, 5g^k M(fc)$k. In Appendix D we show that the perturbative 
formula still works, as long as /nl is not very large. The current observational constraints on /nl already guarantee that 
the perturbative formula for the MFs of the galaxy distribution provides an excellent approximation. 
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C. MFS OF CMB: ANALYTICAL FORMULA VS SIMULATIONS IN THE SACHS-WOLFE LIMIT 

We compare the perturbative formula of MFs for CMB with Monte-Carlo realizations of non-Gaussian temperature 
anisotropy in the Sachs- Wolfe regime, in order to check the accuracy of our formalism. 

The angular power spectrum spectrum is set to be l{l + l)Cf'^ /2n = 10"^°. In the Sachs- Wolfe hmit, AT^^/r = -$/3, 
the non-Gaussian maps of CMB temperature anisotropy may be constructed from the Gaussian maps, ATq/T, by the 
following simple mapping: 



AT^w ATg 

— 7^ — = — 3/nl 



/ATc 



G 




(CI) 



We calculate the MFs from 6000 realizations of the non-Gaussian CMB maps and compare them with the perturbation 
predictions (Eq.[2]). The skewness parameters can be calculated from the reduced bispectrum of CMB in the Sachs- Wolfe 
limit, 

Ohhh = -0/nl(O;^ + + ). (C2) 

Figure C5 shows that the MFs from Monte-Carlo realizations agree with the perturbation predictions very well. The 
comparison with the full simulations that include the full radiation transfer function will be reported elsewhere. (For 
subtlety in this comparison arising from pixclization and boundary effects, see § 4.2.) 

D. VALIDITY OF PERTURBATIVE FORMULAE FOR n— TH ORDER CORRECTIONS OF PRIMORDIAL NON-GAUSSIANITY 

We consider the primordial non-Gaussianity extended to n— th order corrections of a primordial potential field; 

,2f(l) ,3f(2) ,4f(3) ,nAn-l) 

-f = ec^o + 4p(0g - m + 4r0^ + 4r - (0^)) + • ■ ■ + ^^^(^0 - m + ■■■, m 

where e(j)o is an auxiliary random-Gaussian field and represents the coefficients of n-th order of (f>. For convenience, 
we separate the random-Gaussian field into two parts; the e represents the amplitude of the primordial potential power 
spectrum which has the order of 10~* and thereby the fluctuation of (po is the order unity. The parameter /nl in the 
equation (1) corresponds to /nl/2. 



(D3) 
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The Fourier transform of $ is written by 

$ = e^o + 4r^o + ^^^^ 

where $o"''('*) is the Fourier transform of the n-th order term, — (^q)' given by 

^^"^(k) = (2^^^ l^ki • • • y rfk„_i^S(ki) • • • ^S(kn-i)^o(k + ki + • • •k„_i) - (27r)^5c(k)(.^^) 

The polyspectra of $, P^"^ > 2), are defined by 

($(ki) • • • $(k„))e = (27r)3fo(ki + • • •k„)p|"^(ki,. • •,k„) (D4) 

where and correspond to the power spectrum and the bispectrum of $ respectively. 

According to the diagrammatic method by Matsubara (1995), the lowest order terms of e for the connected part of n-th 
order moment, called curnulants, are the n — 1 products of the quadratic moment as follows; 

4'' = ^^InI (Pfo iki)P^o ih) + (eye.)) (D6) 
Pi'^ = eV^\l'(P0o(Ai)^0o(fc2)P0o(|ki+k3|) + (sym.)(12)) 

+ ey^l{P^MPMP4'o{ks) + {cyc.)) (D7) 
Pi'^ = <^'f^T{P't>oiki)P4>oik2)P4>o{\^i + k3|)P^o(|k2 + k4|) + (sym.)(60)) 

+ eV4L/NL(^0o(fci)^0o(fc2)P^o(fc3)P0o(|ki + k4|) + (sym.)(60)) 

+ eV^L iP4>o iWo (fc2)P0o iWo {k4) + (eye.)) (D8) 



p^n) ^ ^2n-2 ^ [ \[ f^' 

di,d2,...,dn \m—l,...,n 



(D9) 



Y[ P^„(|k^ + kB|) + (sym.) 

_kAB(di,d2,....d„,) 

where (sym.)(n) means the addition of n terms with the subscripts symmetric to the previous term and the edge {AB) is 
one of the edges in a tree graph {di,d2, rf„) which satisfies the condition that di + 2d2 + • • • + nd„ = n — 2. 

We obtain the n— th polyspectrum of ^L,k, -Pl"^ t)y 

P^"^(ki,- • •,k„_i) = M{h)M{k2) ■ ■ ■ M(fc„)P|"\ (DIO) 
The n-th order terms of the perturbative formula (eq.[l]) are represented by 

(5E)c (^r'V^5L)e (5r'(V<5L-V5L)V25L)c ^^^^^ 



(^Siyn/2^ (<5E)"/2-i((vJl) . (VJl))' (5D"/'-'((V<5l) • (V<5l))2'- 

These terms are n-th order cumulants of the products of ^l, V^l, and V^(5l divided by n/2 times product of the corre- 
sponding combination of the second moments {Sf) , ((V(5l) • (V(5l)). 

The n-th order cumulants ((5£)c are obtained by the inverse Fourier-transform of Pj^^ as 

{Sl)c = / dkPi'\k)M{krW{kRr (D12) 

{SZ)c = ^2Tr)3n-3 / ^ki • • • J dkn-lPt^M{kr)W{k,R) ■ ■ ■ M{kn)W{knR) (D13) 

The n-th order term of the perturbative formula (eq.[l]) has the following order; 

;<-n\ di+2d-i + -- -ndn=n / \ 

^^^^ -e"-^ ^ n /nl''" ■ (D14) 



\"L/c di,d2,...,dn \m=l,...,n / 

The other n-th order terms in equation (Dll) have the same order of e as (^l)c/(<5l)c^^- 

The above equation is different from the well-known hierarchical condition from the gravitational evolution; 

{Dc-ifr-'- (D15) 

Indeed, the skewness parameters due to the primordial non-Gaussianity have a scale dependence of M{k). Nevertheless, 
the perturbation works well as long as equation (D14) is much smaller than unity, which corresponds to 

|/4"L^|«e-"~10^" (D16) 

Recent observations represented by WMAP (Komatsu et al. 2003) gave constraints on |/nl| < 0(10^). Standard inflation 
models predict that higher-order coefficients are the same order as /nl and thus the perturbation is applicable to the 
actual primordial non-Gaussianity. 



